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Abstract 

THERMUS is a package of C++ classes and functions allowing statistical-thermal 
model analyses of particle production in relativistic heavy-ion collisions to be per- 
formed within the ROOT framework of analysis. Calculations are possible within 
three statistical ensembles; a grand-canonical treatment of the conserved charges 
B, S and Q, a fully canonical treatment of the conserved charges, and a mixed- 
canonical ensemble combining a canonical treatment of strangeness with a grand- 
canonical treatment of baryon number and electric charge. THERMUS allows for 
the assignment of decay chains and detector efficiencies specific to each particle 
yield, which enables sensible fitting of model parameters to experimental data. 
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Nature of problem: 

Statistical-thermal model analyses of heavy-ion collision data require the calculation 
of both primordial particle densities and contributions from resonance decay. A set 
of thermal parameters (the number depending on the particular model imposed) 
and a set of thermalised particles, with their decays specified, is required as input 
to these models. The output is then a complete set of primordial thermal quantities 
for each particle, together with the contributions to the final particle yields from 
resonance decay. 

In many applications of statistical-thermal models it is required to fit experimen- 
tal particle multiplicities or particle ratios. In such analyses, the input is a set of 
experimental yields and ratios, a set of particles comprising the assumed hadron res- 
onance gas formed in the collision and the constraints to be placed on the system. 
The thermal model parameters consistent with the specified constraints leading to 
the best-fit to the experimental data are then output. 

Solution method: 

THERMUS is a package designed for incorporation into the ROOT [2] framework, 
used extensively by the heavy-ion community. As such, it utilises a great deal of 
ROOT's functionality in its operation. ROOT features used in THERMUS include 
its containers, the wrapper TMinuit implementing the MINUIT fitting package, and 
the TMath class of mathematical functions and routines. Arguably the most useful 
feature is the utilisation of CINT as the control language, which allows interactive 
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access to the THERMUS objects. Three distinct statistical ensembles are included 
in THERMUS, while additional options to include quantum statistics, resonance 
width and excluded volume corrections are also available. 

THERMUS provides a default particle list including all mesons (up to the if|(2045)) 
and baryons (up to the listed in the July 2002 Particle Physics Booklet [3]. For 
each typically unstable particle in this list, THERMUS includes a text-file listing 
its decays. With thermal parameters specified, THERMUS calculates primordial 
thermal densities either by performing numerical integrations or else, in the case 
of the Boltzmann approximation without resonance width in the grand-canonical 
ensemble, by evaluating Bessel functions. Particle decay chains are then used to 
evaluate experimental observables (i.e. particle yields following resonance decay). 
Additional detector efficiency factors allow fine-tuning of the model predictions to 
a specific detector arrangement. 

When parameters are required to be constrained, use is made of the 'Numerical 
Recipes in C [1] function which applies the Broyden globally convergent secant 
method of solving nonlinear systems of equations. Since the NRC software is not 
freely-available, it has to be purchased by the user. THERMUS provides the means 
of imposing a large number of constraints on the chosen model (amongst others, 
THERMUS can fix the baryon-to-charge ratio of the system, the strangeness den- 
sity of the system and the primordial energy per hadron). 

Fits to experimental data are accomplished in THERMUS by using the ROOT 
TMinuit class. In its default operation, the standard % 2 function is minimised, 
yielding the set of best-fit thermal parameters. THERMUS allows the assignment 
of separate decay chains to each experimental input. In this way, the model is able 
to match the specific feed-down corrections of a particular data set. 

Running time: Depending on the analysis required, run-times vary from seconds 
(for the evaluation of particle multiplicities given a set of parameters) to several 
minutes (for fits to experimental data subject to constraints). 
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1 Introduction 



The statistical-thermal model has proved extremely successful [T112T3"] in de- 
scribing the hadron multiplicities observed in relativistic collisions of both 
heavy-ions and elementary particles. The methods used in calculating these 
yields have been extensively reviewed in recent years [415] . The success of these 
models has led to the creation of several software codes [6lT7f8] that use ex- 
perimental particle yields as input and calculate the corresponding chemical 
freeze-out temperature (T) and baryon chemical potential (/x.b). In this pa- 
per we present THERMITS, a package of C++ classes and functions, which 
is based on the object-oriented ROOT framework [9J. All THERMITS C++ 
classes inherit from the ROOT base class TObj ect. This allows them to be fully 
integrated into the interactive ROOT environment, allowing all of the ROOT 
functionality in a statistical-thermal model analysis. Recent applications of 
THERMITS include |2^1Ujll|12|13jl4jl5|16jl7jl8|19| . An on-going effort to ex- 
tend the range of applications of THERMITS has led to several publications 
on fluctuations in statistical models [20112 11122] . 

The paper is structured in the following way. In Section 2 an overview is 
presented of the theoretical model on which THERMITS is based. Section 3 
outlines the structure and functionality of the THERMITS code, while Section 
4 explains the installation procedure. 



2 Overview of the Statistical-Thermal Model of Heavy-Ion Colli- 
sions 

2.1 Choice of Ensemble 

Within the statistical-thermal model there is a freedom regarding the ensemble 
with which to treat the quantum numbers B (baryon number), S (strangeness) 
and Q (charge), which are conserved in strong interactions. The introduction of 
chemical potentials for each of these quantum numbers (i.e. a grand-canonical 
description) allows fluctuations about conserved averages. This is a reasonable 
approximation only when the number of particles carrying the quantum num- 
ber concerned is large. In applications of the thermal model to high-energy 
elementary collisions, such as pp, pp and e + e~ collisions [23,24], a canonical 
treatment of each of the quantum numbers is required. Within such a canon- 
ical description, quantum numbers are conserved exactly. In small systems or 
at low temperatures (more specifically, low VT 3 values), a canonical treat- 
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ment leads to a suppression of hadrons carrying non-zero quantum numbers, 
since these particles have to be created in pairs. In heavy-ion collisions, the 
large number of baryons and charged particles generally allows baryon number 
and charge to be treated grand-canonically. However, at the low temperatures 
of the GSI SIS, the resulting low production of strange particles requires a 
canonical treatment of strangeness [25J. This is the so-called mixed-canonical 
approach. 

In order to calculate the thermal properties of a system, one starts with an 
evaluation of its partition function. The form of the partition function obvi- 
ously depends on the choice of ensemble. In the following sections, we consider 
the three ensembles most widely used in applications of the statistical-thermal 
model. 



2.1.1 The Grand- Canonical Ensemble 

This ensemble is the most widely used in applications to heavy-ion colli- 
sions [oT2"6 r f2"7 28 29 30 31 32 33 l3"Hl3l)] . Within this ensemble, conservation laws 
for energy and quantum or particle numbers are enforced on average through 
the temperature and chemical potentials. 

In the case of a multi-component hadron gas of volume V and temperature 
T, the logarithm of the total partition function is given by, 



\nZ GC (T,V,{^}) = E £^ d 3 p\n(l±e-^-^f\ (1) 



where and are, respectively, the degeneracy and chemical potential of 
hadron species i, f3 = 1/T, while = ^Jp 2 + m 2 , where m; is the particle 
mass. The plus sign refers to fermions and the minus sign to bosons. 

Since in relativistic heavy-ion collisions it is not individual particle numbers 
that are conserved, but rather the quantum numbers B, S and Q, the chemical 
potential for particle species i is given by, 

Vi = B~i[iB + Sifis + QifiQ, (2) 

where 5,, Si and Qi are the baryon number, strangeness and charge, respec- 
tively, of hadron species i, and /is and /ig are the corresponding chemical 
potentials for these conserved quantum numbers. 



species i 
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Once the partition function is known, the particle multiplicities, entropy and 
pressure are obtained by differentiation: 



N? C = T— , 3 

S GC = -^(TlnZ GC ), (4) 
rr d\nZ GG 

p = T ^T- ^ 

Furthermore, the energy is given by, 



E GC = T*^-+ J2 »iNf C - (6) 

species i 

Using the prescription for the particle multiplicity, 



N?° = t(Tl) k+1 ^K 2 (^j e?*» = ± W", (7) 

where we have introduced the zf. Similar expressions exist for the energy, en- 
tropy and pressure. 



In practice, the Boltzmann approximation (i.e. retaining just the k = 1 term 
in Equation ((Tj)) is reasonable for all particles except the pions. In this ap- 
proximation, 



\nZ GG (T,V,{^}) = E J d 3 p e~^-^ = £ z\ eP", (8) 

species i V J species i 

where z\ is the single-particle partition function of hadron species i. Further- 
more, under this approximation, P = ^species i N GG T/V both for massive and 
massless particles, which is certainly not true for quantum statistics. 



Since the use of quantum statistics requires numerical integration (or evalu- 
ation of infinite sums), while Boltzmann statistics can be implemented ana- 
lytically, it is worthwhile to identify those regions in which quantum statis- 
tics deviate greatly from Boltzmann statistics. In most applications of the 
statistical-thermal model, only a small region of the \i — T parameter space is 
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of interest. Using the freeze-out condition of constant E/N [36], the thermal 

parameters, and hence the percentage deviation from Boltzmann statistics, 
can be determined as a function of the collision energy yfs [37] . From such an 
analysis it is evident that, for pions, quantum statistics must be implemented 
at all but the lowest energies (deviation at the level of 10%), while, for kaons, 
the deviation peaks at between 1 and 2%. For all other mesons, the deviation 
is below the 1% level. For baryons, the deviation is extremely small for all 
except the protons at small y/s. 

When quantum statistics are applied, restrictions have to be imposed on the 
chemical potentials so as to avoid Bose-Einstein condensation. The Bose- 
Einstein distribution function diverges if, 



Such Bose-Einstein condensation is avoided, provided that the chemical po- 
tentials of all bosons included in the resonance gas are less than their masses 
(i.e. (it < rrii). 



2.1.2 The Canonical Ensemble 

Within this ensemble, quantum number conservation is exactly enforced. Con- 
sidering the fully canonical treatment of B, S and Q in the Boltzmann ap- 
proximation, as investigated in [35], the partition function for the system is 
given by, 



(9) 




— 7T — 7T 



x exp 2 



mesons j 



E 




(10) 



where, 



E 



1 e i (Sj ( l>s+Qj<i>Q) 



baryons j 
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Zq represents the contribution of those hadrons with no net charges, and the 
sums over mesons and baryons extend only over the particles (i.e. not the 
anti-particles). 



Once the partition function is known, we can calculate all thermodynamic 
properties of the system. Using thermodynamic relations it follows that, 



S = -^(TlnZ BAQ ), 



and, 



P = T 



d In Z BtS ,Q 
dV ' 



(12) 



Furthermore, the multiplicity of hadron species i within this ensemble, N^' s,< ^ , 
is calculated by multiplying the single-particle partition function for particle 
i, appearing in the canonical partition function, by a fictitious fugacity Aj, 
differentiating with respect to A,, and then setting A, to 1: 



b,s,q = d\nZ BjS ,Q(\i 



Ai=l 



(13) 



Following these prescriptions, 



N, 



B,S,Q 



'B-Bi,S-Si,Q-Qi 
Z B,S,Q 



N. 



GC 



(14) 



S B ^ = \nZ BAQ + £ ( 

species i \ 

pB,S,Q = ^ ( 
species i \ 

E B,S,Q = ( 
species i \ 



Z R_ 



B—Bi,S—Si,Q—C 



J B,S,Q 



species % \ B,S,Q 

z 



B-Bi,S-Si,Q-Qi \ pGC 



E? C 



T 



B- Bi,S-Si,Q-Qi pGC 



M=0 ' 



(15) 
(16) 
(17) 



8 



One notices that, in the Boltzmann approximation, the particle and energy 
density and pressure of particle species i, within the canonical ensemble, differ 
from that in the grand-canonical formalism, with all chemical potentials set to 
zero, by a multiplicative factor (Zb-b^s-s^q-qJ 'Zb,s,q)- This correction fac- 
tor depends only on the thermal parameters of the system and the quantum 
numbers of the particle (i.e. the correction for the A + and p are the same). 
The entropy is, however, slightly different; the total entropy cannot be split 
into the sum of contributions from separate particles. 



Now, 



lim 



J B-B t ,S~S^Q-Q, 

Zb,s,q 



qB^b/T e Siix s /T e QnJ, Q /T _ 



Thus, for large systems, the grand-canonical results for the particle number, 
entropy, pressure and energy are approached [38J. 



2.1.3 The Mixed- Canonical (Strangeness-Canonical) Ensemble 

Within this ensemble, the strangeness in the system is fixed exactly by its 
initial value of S, while the baryon and charge content are treated grand- 
canonically. For a Boltzmann hadron gas of strangeness S, 



Zft = — 



2tt 

— 7T 

x exp 



>s e 



-iS<t> s 



j2 I d* P v. 



-hadrons i 

where the sum over hadrons includes both particles and anti-particles and 



(19) 



(20) 



Applying the same prescription for the evaluation of the particle multiplicities 
as discussed for the canonical ensemble, it follows that, 



Nf 



j s-s t 



N 



GC 



(21) 



Furthermore, 
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S s = \nZ s + J2 

species i 

z 



S-Si 



E 



GC 



-IH N, 



GC 



E 



species i 



species i 



S—Si \ pGC 



MS=0 



(22) 

(23) 
(24) 



As in the case of the canonical ensemble, the strangeness-canonical results, 
in the Boltzmann approximation, differ from those of the grand-canonical 
ensemble, with fi$ = 0, by multiplicative correction factors which depend, in 
this case, only on the thermal parameters and the strangeness of the particle 
concerned. For large systems and high temperatures, these correction factors 
approach the grand-canonical fugacities, i.e., 



lim 



j s-s t 



(25) 



The expression for Z s can be reduced [39] to, 



+oo +CO 



Z S = Z X J2 I \3m+2n-S\(x 1 ) I\ n \(x 2 ) i| m |(x 3 ) 



(26) 



m=— oo n=— oo 



where Z is the contribution to the total partition function of the non-strange 
hadrons, while, 



Xi = 2 \ k+ik 



and, 



l k +i 
k—i 



with, 



(z = 1,2,3), 



(27) 



(z = 1,2,3), 



(2c 



hadrons j with Sj =m 



V. 



(29) 



In [4"0|4"T] it is suggested that two volume parameters be used within canoni- 
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cal ensembles; the fireball volume at freeze-out, Vf, which provides the overall 
normalisation factor fixing the particle multiplicities from the correspond- 
ing densities, and the correlation volume, V c , within which particles fulfill 
the requirement of local conservation of quantum numbers. In this way, by 
taking V c < Vf, it is possible to boost the strangeness suppression. In fact, 
this was shown to be required to reproduce experimental heavy-ion collision 
data [30121] . 



2.2 Feeding from Unstable Particles 

Since the particle yields measured by the detectors in collision experiments in- 
clude feed-down from heavier hadrons and hadronic resonances, the primordial 
hadrons are allowed to decay to particles considered stable by the experiment 
before model predictions are compared with experimental data. For example, 
the total 7r + yield is given by, 

N, + = £ N^Br{i^^), (30) 

species i 

where Br(i — > tt + ) is the number of 7r + 's into which a single particle of species 
i decays. As shown in Figure dj approximately 70% of 7r + 's originate from 
resonance decay at RHIC energies. Thus, a full treatment of resonances is es- 
sential in any statistical-thermal analysis. 

The inclusion of a mass cut-off in the measured resonance mass spectrum is 
motivated by the realisation that the time scale of a relativistic collision does 
not allow the heavier resonances to reach chemical equilibrium [12]. This as- 
sumes that inelastic collisions drive the system to chemical equilibrium. If the 
hadronisation process follows a statistical rule, then all resonances should, in 
principle, be included [3Tj . This is problematic, since data on the heavy res- 
onances is sketchy. The situation is saved by the finite energy density of the 
system, resulting in a chemical freeze-out temperature at RHIC of approxi- 
mately 160-170 MeV [T3l 1 1 . which strongly suppresses these heavy resonances 
and justifies their exclusion from the model. It is, however, important to check 
the sensitivity of the extracted thermal parameters to the chosen cut-off. The 
effect of the mass cut-off on particle ratios was studied in [37] . 

The finite width of the resonances is especially important at the low tempera- 
tures of the SIS. Resonance widths are included in the thermal model by dis- 
tributing the resonance masses according to Breit-Wigner forms [23i24ll27ll25f42 
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Fig. 1. The energy dependence of the decay contribution to the final n + yield as pre- 
dicted by the statistical-thermal model. Calculations performed within the grand- 
canonical formalism, assuming the constant E/N freeze-out criterion |36j (weak 
decays excluded). 

This amounts to the following modification in the integration of the Boltzmann 
factor 1251: 



J d 3 p exp 



\Jp 2 + m 2 



T 



d 3 p J ds exp 



Vp 2 



T 



7T 



mT 



is — m 2 ) 2 + m 2 T 2 ' 



(31) 



where Y is the width of the resonance concerned, with threshold limit m t h res hoid 
and mass m, and y/s is integrated over the interval [m - 5m, m + 2T], where 
5m = min[m - m thresho , d , 2T]. 



2. 3 Deviations from Equilibrium Levels 



The statistical-thermal model applied to elementary e + e~, pp and pp colli- 
sions [2"3"1T2"^] indicates the need for an additional parameter, 75 (first intro- 
duced as a purely phenomenological parameter pESP?] ). to account for the 
observed deviation from chemical equilibrium in the strange sector. Since a 
canonical ensemble was considered in these analyses, there is an additional 
strangeness suppression at work, on top of the canonical suppression. Although 
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strangeness production is expected to be greatly increased in AA collisions, 
due to the larger interaction region and increased hadron rescattering, a num- 
ber of recent analyses [45,48,4 91l5UII51f52| have found such a factor necessary 
to accomplish a satisfactory description of data. 

Allowance for possibly incomplete strangeness equilibration is made by multi- 
plying the Boltzmann factors of each particle species in the partition function 

I S ' I 

(or thermal distribution function fi(x,p)) by t s , where \Si\ is the number 
of valence strange quarks and anti-quarks in species i (for example, for the 
meson, with an ss pair, |S^| = 2). The value 75 = 1 obviously corresponds 
to complete strangeness equilibration. 

It has been suggested [53] that a similar parameter, j q , should be included in 
thermal analyses to allow for deviations from equilibrium levels in the non- 
strange sector. Furthermore, as collider energies increase, so does the need for 
the inclusion of charmed particles in the statistical-thermal model, with their 
occupation of phase-space possibly governed by an additional parameter, jc- 



2.4 Excluded Volume Corrections (Grand- Canonical Ensemble) 

At very high energies, the ideal gas assumption is inadequate. In fact, the total 
particle densities predicted by the thermal model, with parameters extracted 
from fits to experimental data, far exceed reasonable estimates and measure- 
ments based on yields and the system size inferred by pion interferometry |54j . 
It becomes necessary to take into account the Van der Waals-type excluded 
volume procedure [5"4"|55|l56j . At the same fixed T and /xg, all thermodynamic 
functions of the hadron g clS clXG smaller than in the ideal hadron gas, and 
strongly decrease with increasing excluded volume. 

Van der Waals-type corrections are included by making the following substi- 
tution for the volume V in the canonical (with respect to particle number) 
partition function, 

V^V- "i N i> ( 32 ) 

hadrons i 

where Ni is the number of hadron species i, and Vi = 4 (4/37rrf ) is its proper 
volume, with r.j its hard-sphere radius. This then leads to the following tran- 
scendental equation for the pressure of the gas in the grand-canonical ensemble 
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(assuming h particle species): 



P( T , --■> t J 'h) = J2 P i deal ( T > &)> (33) 
i=i 

with, 

fii = Hi-v i P(T > fjii,...,Hh)- (34) 
The particle, entropy and energy densities are given by, 



TH{T,p 1 ,...,p h ) = —— ideal , 35 

1 + Ej ^n* (T, /ij) 

^^'■^^l-fE,^"^)' (36) 



and, 



1 + E J -^nf eai (T,/x j ) : 



e{T,m,...,n h ) = - ■ - - — , (37) 



respectively One sees that two suppression factors enter. The first suppression 
is due to the shift in chemical potential /ij — ► /Ij. In the Boltzmann approxima- 
tion, this leads to a suppression factor e~ ViP l T in all thermodynamic quantities. 
The second suppression is due to the [1 + Ej P^ al {T, fij)}' 1 factor. 



In ratios of particle numbers, although the denominator correction cancels 
out, the shift in chemical potentials leads to a change in the case of quantum 
statistics. In the Boltzmann case, even these corrections cancel out, provided 
that the same proper volume parameter v is applied to all species. 



3 The Structure of THERMUS 



3. 1 Introduction 



The three distinct ensemble choices outlined in Sections I2.1.1tl27l~3l are im- 
plemented in THERMUS. As input to the various thermal model formalisms 
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one needs first a set of particles to be considered thermalised. When combined 
with a set of thermal parameters, all primordial densities (i.e. number density 
as well as energy and entropy density and pressure) are calculable. Once the 
particle decays are known, sensible comparisons can be made with experimen- 
tally measured yields. 



In THERMUS, the following units are used for the parameters: 



Parameter Unit 

Temperature (T) GeV 

Chemical Potential (/i) GeV 

Radius fm 



Quantities frequently output by THERMUS are in the following units: 



Quantity Unit 
Number Densities (n) fm -3 

Energy Density (e) GeV.fm -3 
Entropy Density (s) fm -3 

Pressure (P) GeV.fm -3 
Volume (V) fm 3 



In the subsections to follow, we explain the basic structure and functionality 
of THERMUS (shown diagrammatically in Figure [2]) by introducing the ma- 
jor THERMUS classes in a bottom-up approach. We begin with a look at the 
TTMParticle object {3 



2 It is a requirement that all ROOT classnames begin with a 'T'. THERMUS 
classnames begin with 'TTM' for easy identification. 
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Particle 

TTMParticle 




Parameter Set 

TTMParameterSet 



Particle properties 
Particle decays 



Thermal Particle 

TTMThermalParticle 



A 



Functions to calculate primordial 
particle-, energy- and entropy densities 
and pressure 



Particle Set 

TTMParticleSet 



Functions and 
structures used in 
thermal model classes 



Thermal Model 

TTMThermalModel 




Primordial properties and 
feed-down contributions of 
all constituents calculable 
Total fireball properties 
now known 



Densities 

TTMDensOb j 



Experimental Data 

TTMYield 



Container for thermal 
densities output by 
thermal model classes 



Thermal Fit 

TTMThermalFit 



Best-fit thermal parameters 
can now be determined 



Fig. 2. The basic structure of THERMUS (only the most fundamental base classes 
are shown) . 



3.2 The TTMParticle Class 



The properties of a particle applicable to the statistical-thermal model are 
grouped in the basic TTMParticle object: 



********* LISTING FOR PARTICLE Delta(1600)0 ********* 

ID = 32114 
Deg. = 4 
STAT = 1 

Mass = 1.6 GeV 

Width = 0.35 GeV 

Threshold = 1.07454 GeV 

Hard sphere radius = 

B = 1 
S = 
Q = 
Charm = 
Beauty = 
Top = 

UNSTABLE 

Decay Channels: 

Summary of Decays : 



Besides the particle name, 'Delta(1600)0' in this case, its Particle Data Group 
(PDG) numerical ID is also stored. This provides a far more convenient means 
of referencing the particle. The particle's decay status is also noted. In this 
case, the A(1600)° is considered unstable. Particle properties are input using 
the appropriate setters. 
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3.2.1 Inputting and Accessing Particle Decays 

The TTMParticle class allows also for the storage of a particle's decays. These 
can be entered from file. As an example, consider the decay file of the A(1600)°: 



11.67 


2112 


111 


5.83 


2212 


-211 


29.33 


2214 


-211 


3.67 


2114 


111 


22. 


1114 


211 


8.33 


2112 


113 


4.17 


2212 


-213 


15. 


12112 


111 


7.5 


12212 


-211 



Each line in the decay file corresponds to a decay channel. The first column 
lists the branching ratio of the channel, while the subsequent tab-separated 
integers represent the PDG ID's of the daughters (each line (channel) can 
contain any number of daughters). The decay channel list of a TTMParticle 
object is populated with TTMDecayChannel objects by the SetDecayChannels 
function, with the decay file the first argument (only that part of the output 
that differs from the previous listing of the particle information is shown) 
(Note: in the example below $THERMUS must be entered in full): 

root [ ] part->SetDecayChannels("$THERMUS/particles/Delta\(1600\) 0_decay.txt") 
root [ ] part->List () 



********* LISTING FOR PARTICLE Delta(1600)0 ********* 



UNSTABLE 



Decay Channels : 
BRatio: 0.1167 
BRatio: 0.0583 
BRatio: 0.2933 
BRatio: 0.0367 
BRatio: 0.22 



Daughters : 
Daughters : 
Daughters : 
Daughters : 
Daughters : 
Daughters : 
Daughters : 



2212 
2214 
2114 
1114 



2112 



111 

-211 

-211 



211 



111 



BRatio: 0.0833 
BRatio: 0.0417 



2112 
2212 



113 
-213 
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BRatio: 0.15 
BRatio: 0.075 



Daughters: 12112 111 

Daughters: 12212 -211 



O LLllllllCLX V Ul 


\* CL V O * 






in 


30 . 34% 


2212 


10% 


-211 


42 . 66% 


2214 


29 . 33% 


2114 


3 . 67% 


1114 


22% 


211 


22% 


113 


8 . 33% 


-213 


4 . 17% 


12112 


15% 


12212 


7.5% 



In many cases, the branching ratios of unstable hadrons do not sum to 100%. 
This can, however, be enforced by scaling all branching ratios. This is achieved 
when the second argument of SetDecayChannels is set to true (it is false by 
default). 

In addition to the list of decay channels, a summary list of TTMDecay objects 
is generated in which each daughter appears only once, together with its total 
decay fraction. This summary list is automatically generated from the decay 
channel list when the SetDecayChannels function is called. As an example, 
the summary list of the A + contains the following entries: p: 2/3, n: 1/3, 7r + : 
1/3, vr°: 2/3. 

An existing TList can be set as the decay channel list of the particle, using 
the SetDecayChannels function. This function calls UpdateDecaySummary, 
thereby automatically ensuring consistency between the decay channel and 
decay summary lists. 

The function SetDecayChannelEf f iciency sets the reconstruction efficiency 
of the specified decay channel to the specified percentage (it has a default 
value of 100%). Again, a consistent decay summary list is generated. 
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Access to the TTMDecayChannel objects in the decay channel list is achieved 
through the TTMDecayChannel* GetDecayChannel method. If the extracted 
decay channel is subsequently altered, UpdateDecaySummary must be called 
to ensure consistency of the summary list. 



3.3 The TTMParticleSet Class 

The thermalised fireballs considered in statistical-thermal models typically 
contain approximately 350 different hadron and hadronic resonance species. 
To facilitate fast retrieval of particle properties, the TTMParticle objects of 
all constituents are stored in a hash table in a TTMParticleSet object. Other 
data members of this TTMParticleSet class include the filename used to in- 
stantiate the object and the number of particle species. Access to the entries 
in the hash table is through the PDG ID's. 



3.3.1 Instantiating a TTMParticleSet Object 

In addition to the default constructor, the following constructors exist: 

TTMParticleSet *set = new TTMParticleSet (char *f ile) ; 
TTMParticleSet *set = new TTMParticleSet (TDatabasePDG *pdg) ; 

The first constructor instantiates a TTMParticleSet object and inputs the 
particle properties contained in the specified text file. As an example of such a 
file, $THERMUS/particles/PartList_PPB2002 . txt contains a list of all mesons 
(up to the ^(2045)) and baryons (up to the Or) listed in the July 2002 Par- 
ticle Physics Booklet [57] (195 entries). Only particles need be included, since 
the anti-particle properties are directly related to those of the corresponding 
particle. The required file format is as follows: 

Delta(1600)0 32114 4 +1 1.60000 1 

0.35000 1.07454 (npiO) 

• stability flag (1 for stable, for unstable) 

• particle name 

• PDG ID (used for all referencing) 

• spin degeneracy 

• statistics (+1 for Fermi-Dirac, -1 for Bose-Einstein, for Boltzmann) 
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Fig. 3. The mass distribution of the resonances included in PartList_PPB2002.txt. 



• mass in GeV 

• strangeness 

• baryon number 

• charge 

• absolute strangeness content \S\ — #s + #s (e.g., \S^\ = 2) 

• width in GeV 

• threshold in GeV 

• string recording the decay channel from which the threshold is calculated if 
the particle's width is non-zero 

All further particle properties have to be set with the relevant setters (e.g. 
the charm, absolute charm content and hard-sphere radius). By default, all 
properties not listed in the particle list file are assumed to be zero. 



Figure [3] shows the distribution of resonances (both particle and anti-particle) 
derived from $THERMUS/particles/PartList_PPB2002.txt. As collider en- 
ergies increase, so does the need to include also the higher mass resonances. 
Although the TTMParticle class allows for the properties of charmed parti- 
cles, these particles are not included in the default THERMUS particle list. 
If required, these particles have to be input by the user. The same applies to 
the hadrons composed of b and t quarks. 
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It is also possible to use a TDatabasePDG object to instantiate a particle se 
TDatabasePDG objects also read in particle information from text files. The 
default file is $ROOTSYS/etc/pdg_table .txt and is based on the parameters 
used in PYTHIA ^E\. 

The constructor TTMParticleSet (TDatabasePDG *pdg) extracts only those 
particles in the specified TDatabasePDG object in particle classes 'Meson', 
'CharmedMeson', 'HiddenCharmMeson', 'B-Meson', 'Baryon', 'CharmedBaryon' 
and 'B-Baryon', as specified in $ROOTSYS/etc/pdg_table .txt, and includes 
them in the hadron set. Anti-particles must be included in the TDatabasePDG 
object, as they are not automatically generated in this constructor of the 
TTMParticleSet class. 

The default file read into the TDatabasePDG object, however, is incomplete; the 
charm, degeneracy, threshold, strangeness, 15*1, beauty and topness of the par- 
ticle are not included. Although the TDatabasePDG: :ReadPDGTable function 
and default file allow for isospin, J3, spin, flavor and tracking code to be entered 
too, the default file does not contain these values. Furthermore, all particles 
are made stable by default. Therefore, at present, using the TDatabasePDG 
class to instantiate a TTMParticleSet class should be avoided, at least until 
pdg_table. txt is improved. 

3.3.2 Inputting Decays 

Once a particle set has been defined, the decays to the stable particles in 
the set can be determined. After instantiating a TTMParticleSet object and 
settling on its stable constituents (the list of stable particles can be modified 
by adjusting the stability flags of the TTMParticle objects included in the 
TTMParticleSet object), decays can be input using the InputDecays method. 
Running this function populates the decay lists of all unstable particles in the 
set, using the decay files listed in the directory specified as the first argument. 
If a file is not found, then the corresponding particle is set to stable. For each 
typically unstable particle in $THERMUS/particles/PartList_PPB2002.txt, 
there exists a file in $THERMUS/particles listing its decays. The filename 
is derived from the particle's name (e.g. Delta(1600) O_decay.txt for the 
A(1600)°). There are presently 195 such files, with entries based on the Par- 
ticle Physics Booklet of July 2002 [57]. The decays of the corresponding 
anti-particles are automatically generated, while a private recursive function, 
GenerateBRatios, is invoked to ensure that only stable particles feature in 
the decay summary lists. The second argument of InputDecays, when set to 

3 In order to have access to TDatabasePDG and related classes, one must first load 
$RDDTSYS/lib/libEG . so 
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true, scales the branching ratios so that their sum is 100%. As an example, 
consider the following (again only part of the listing is shown) (Note: in the 
example below $THERMUS must be entered in full): 



root [ ] TTMParticleSet set("$THERMUS/particles/PartList_PPB2002.txt") 

root [ ] set . InputDecays("$THERMUS/particles/" .true) 

root [ ] TTMParticle *part = set . GetParticle (32114) 

root [ ] part->List () 

********* LISTING FOR PARTICLE Delta(1600)0 ********* 



UNSTABLE 



Decay Channels: 



BRatio: 


0. 


.108558 


Daughters : 


2112 


111 


BRatio: 


0. 


,0542326 


Daughters : 


2212 


-211 


BRatio: 


0. 


.272837 


Daughters : 


2214 


-211 


BRatio: 


0. 


.0341395 


Daughters : 


2114 


111 


BRatio: 


0. 


.204651 


Daughters : 


1114 


211 


BRatio: 


0. 


. 0774884 


Daughters : 


2112 


113 


BRatio: 


0. 


.0387907 


Daughters : 


2212 


-213 


BRatio: 


0. 


.139535 


Daughters : 


12112 


111 


BRatio: 


0. 


.0697674 


Daughters : 


12212 


-211 



Summary of Decays : 



2112 60.6774% 

111 62.5704% 

2212 39.3226% 

-211 83.9999% 

211 44.6773% 



For particle sets based on TDatabasePDG objects, decay lists should be pop- 
ulated through the function InputDecays (TDatabasePDG *). This function, 
however, does not automatically generate the anti-particle decays from those 
of the particle. Instead, the anti-particle decay list is used. Since the decay 
list may include electromagnetic and weak decays to particles other than the 
hadrons stored in the TTMParticleSet object, each channel is first checked 
to ensure that it contains only particles listed in the set. If not, the channel 
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is excluded from the hadron's decay list used by THERMUS. As mentioned 
earlier, care should be taken when using TDatabasePDG objects based on the 
default file, as it is incomplete. 

An extremely useful function is ListParents (Int_t id), which lists all of 
the parents of the particle with the specified PDG ID. This function uses 
GetParents (TList *parents, Int_t id), which populates the list passed 
with the decays to particle id. Note that these parents are not necessarily 
'direct parents'; the decays may involve unstable intermediates. 

3.3.3 Customising the Set 

The AddParticle and RemoveParticle functions allow customisation of par- 
ticle sets. Particle and anti-particle are treated symmetrically in the case of 
the former; if a particle is added, then its corresponding anti-particle is also 
added. This is not the case for the RemoveParticle function, however, where 
particle and anti-particle have to be removed separately. 

Mass-cuts can be performed using MassCut (Double_t x) to exclude all hadrons 
with masses greater than the argument (expressed in GeV). Decays then have 
to be re-inserted, to remove the influence of the newly-excluded hadrons from 
the decay lists. 

The function SetDecayEf f iciency allows the reconstruction efficiency of the 
decays from a specified parent to the specified daughter to be set. Changes 
are reflected only in the decay summary list of the parent (i.e. not the decay 
channel list). Note that running UpdateDecaySummary or GenerateBRatios 
will remove any such changes, by creating again a summary list consistent 
with the channel list. 

In addition to these operations, users can input their own particle sets by 
compiling their own particle lists and decay files. 



3.4 The TTMParameter Class 

This class groups all relevant information for parameters in the statistical- 
thermal model. Data members include: 
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f Name 


- the parameter name, 


i va±ue 


- the parameter value, 


f Error 


- the parameter error, 


fFlag 


- a flag signalling the type of parameter (constrain, fit, 




fixed, or uninitialised), 


f Status 


- a string reflecting the intended treatment or action taken 



In addition to these data members, the following, relevant to fit-type param- 
eters, are also included: 

f Start - the starting value in a fit, 

f Min - the lower bound of the fit-range, 

f Max - the upper bound of the fit-range, 

f Step - the step-size. 

The constructor and SetParameter (TString name, Double_t value, Double_t 
error) function set the parameter to fixed-type, by default. The parameter- 
type can be modified using the Constrain, Fit or Fix methods. 



3.5 The TTMParameterSet Class 

The TTMParameterSet class is the base class for all thermal parameter set 
classes. Each derived class contains its own TTMParameter array, with size 
determined by the requirements of the ensemble. The base class contains a 
pointer to the first element of this array. In addition, it stores the constraint 
information. 

All derived classes contain the function double GetRadius. In this way, TTMParameterSet 
is able to define a function, double GetVolume, which returns the volume re- 
quired to convert densities into total fireball quantities. 

TTMParameterSetBSQ, TTMParameterSetBQ and TTMParameterSetCanBSQ are 

the derived classes. 
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3. 5. 1 TTMParame t erSe tBSQ 

This derived class, applicable to the grand-canonical ensemble, contains the 
parameters: 

T fi B fi S fi Q fi C Is lc R, 

where R is the fireball radius, assuming a spherical fireball (i.e. V = 4/3irR 3 ). 
In addition, the B/2Q ratio and charm and strangeness density of the system 
are stored here. In the constructor, all errors are defaulted to zero, as is R, 
fic, S/V, C/V and B/2Q, while 7^ is defaulted to unity. 

Each parameter has a getter (e.g. TTMParameter* GetTPar), which returns 
a pointer to the requested TTMParameter object. In this class, fi$ an d HQ 
can be set to constrain-type using ConstrainMuS and ConstrainMuQ, where 
the arguments are the required strangeness density and B/2Q ratio, respec- 
tively. No such function exists for fie, since constraining functions are not 
yet implemented for the charm density. Each parameter of this class can be 
set to fit-type, using functions such as FitT (where the fit parameters have 
reasonable default values), or fixed-type, using functions such as FixMuB. 

3.5.2 TTMParameterSetBQ 

This derived class, applicable to the strangeness-canonical ensemble (strangeness 
exactly conserved and B and Q treated grand-canonically), has the parame- 
ters: 

T [x B h q 75 R c R, 

where R c is the canonical or correlation radius; the radius inside which strangeness 
is exactly conserved. The fireball radius R, on the other hand, is used to con- 
vert densities into total fireball quantities. In addition, the required B/2Q 
ratio is also stored, as well as the strangeness required inside the correlation 
volume (which must be an integer). 

In addition to the same getters and setters as the previous derived class, it is 
possible to set /iQ to constrain-type by specifying the B/2Q ratio in the argu- 
ment of ConstrainMuQ. The strangeness required inside the canonical volume 
is set through the Sets method. This value is defaulted to zero. The function 
ConserveSGlobally fixes the canonical radius, R c , to the fireball radius, R. 
As in the case of the TTMParameterSetBSQ class, there also exist functions to 
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set each parameter to fit or fixed- type. 



3. 5. 3 TTMParame t erSe t CanBSQ 



This set, applicable to the canonical ensemble with exact conservation of B, 
S and Q, contains the parameters: 



T B S Q 7s R. 



Since all conservation is exact, there are no chemical potentials to satisfy 
constraints. Again, the same getters, setters and functions to set each pa- 
rameter to fit or fixed-type exist, as in the case of the previously discussed 
TTMParameterSet derived classes. 



3.5.4 Example 

As an example, let us define a TTMParameterSetBQ object. By default, all 
parameters are initially of fixed-type. Suppose we wish to fit T and /ib, and 
use fiQ to constrain the B/2Q ratio in the model to that in Pb+Pb collisions: 

root [] TTMParameterSetBQ parBQ (0. 160,0. 2,-0. 01 ,0. 8,6. ,6. ) 

root [ ] parBQ.FitT(0. 160) 

root [ ] parBQ . FitMuB (0.2) 

root [ ] parBQ . ConstrainMuQ ( 1.2683) 

root [ ] parBQ. List () 

***************************** Thermal Parameters **************************** 



Strangeness inside Canonical Volume = 



T 



0.16 



(to be FITTED) 
start : 0.16 
range: 0.05 — 0.18 
step: 0.001 



muB 



0.2 



(to be FITTED) 
start: 0.2 
range : — 
step: 0.001 



0.5 



muQ 



-0.01 



(to be CONSTRAINED) 
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B/2Q: 1.2683 



gammas 



0.8 



(FIXED) 



Can. radius 



6 



(FIXED) 



radius 



6 



(FIXED) 



Parameters unconstrained 



Note the default parameters for the T and \1b fits. Obviously, no constraining 
or fitting can take place yet; we have simply signalled our intent to take these 
actions at some later stage. 

3.6 The TTMThermalParticle Class 

By combining a TTMParticle and TTMParameterSet object, a thermal particle 
can be created. The TTMThermalParticle class is the base class from which 
thermal particle classes relevant to the three currently implemented thermal 
model formalisms, TTMThermalParticleBSQ, TTMThermalParticleBQ and 
TTMThermalParticleCanBSQ, are derived. Since no particle set is specified, the 
total fireball properties cannot be determined. Thus, in the grand-canonical 
approach, the constraints cannot yet be imposed to determine the values of 
the chemical potentials of constrain- type, while, in the strangeness-canonical 
and canonical formalisms, the canonical correction factors cannot yet be calcu- 
lated. Instead, at this stage, the chemical potentials and/or correction factors 
must be specified. 

Use is made of the fact that, in the Boltzmann approximation, E/V, N/V 
and P, in the canonical and strangeness-canonical ensembles, are simply the 
grand-canonical values, with the chemical potential(s) corresponding to the 
canonically-treated quantum number(s) set to zero, multiplied by a particle- 
specific correction factor. This allows the functions for calculating E/V, N/V 
and P in the Boltzmann approximation to be included in the base class, which 
then also contains the correction factor as a data member (by definition, this 
correction factor is 1 in the grand-canonical ensemble). 

Both functions including and excluding resonance width, T, are implemented 
(e.g. double Dens ityBoltzmannNoWidth and double EnergyBoltzmannWidth) 
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When width is included, a Breit-Wigner distribution is integrated over between 
the limits [max(m — 2T, m thrcshold ), m + 2r]. 



3.6.1 TTMThermalParticleBSQ 

This class is relevant to the grand-canonical treatment of B, S and Q. In 
addition to the functions for calculating E/V, N/V and P in the Boltz- 
mann approximation, defined in the base class, functions implementing quan- 
tum statistics for these quantities exist in this derived class (e.g. double 
EnergyQStatNoWidth and double PressureQStatWidth). Additional mem- 
ber functions of this class calculate the entropy using either Boltzmann or 
quantum statistics, with or without width. 

In the functions calculating the thermal quantities assuming quantum statis- 
tics, it is first checked that the integrals converge for the bosons (i.e. there is no 
Bose-Einstein condensation). The check is performed by the bool ParametersAllowed 
method. A warning is issued if there are problems and zero is returned. 

This class also accommodates charm, since the associated parameter set in- 
cludes He and 7c, while the associated particle may have non-zero charm. 



3.6.2 TTMThermalParticleBQ 

This class is relevant to the strangeness-canonical ensemble. At present, this 
class is only applied in the Boltzmann approximation. Under this assumption, 
N/V, E/V and P are given by the grand-canonical result, with fx s set to zero, 
up to a multiplicative correction factor. Since the total entropy does not split 
into the sum of particle entropies, no entropy calculation is made in this class. 



3.6.3 TTMThermalParticleCanBSQ 

This class is relevant to the fully canonical treatment of B, S and Q. At 
present, as in the case of TTMThermalParticleBQ, this class is only applied 
in the Boltzmann approximation. Also, since the total entropy again does not 
split into the sum of particle entropies, no entropy calculation is made here. 
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3.6.4 Example 



Let us construct a thermal particle, within the strangeness-canonical ensemble, 
from the A(1600)° and the parameter set previously defined. Since this particle 
has zero strangeness, a correction factor of 1 is passed as the third argument 
of the constructor: 

root [ ] TTMThermalParticleBQ therm_delta(part ,&parBQ, 1 . ) 

root [ ] therm_delta.DensityBoltzmannNoWidth() 

(Double_t)8.15072671710089913e-04 

root [ ] therm_delta.EnergyBoltzmannWidth() 

(Double_t)2.29185316377137748e-03 



3.7 The TTMThermalModel Class 

Once a parameter and particle set have been specified, these can be com- 
bined into a thermal model. TTMThermalModel is the base class from which 
the TTMThermalModelBSQ, TTMThermalModelBQ and TTMThermalModelCanBSQ 
classes are derived. A string descriptor is included as a data member of the 
base class to identify the type of model. This is used, for example, to handle 
the fact that the number of parameters in the associated parameter sets is 
different, depending on the model type. 

All derived classes define functions to calculate the primordial particle, energy 
and entropy densities, as well as the pressure. These thermal quantities are 
stored in a hash table of TTMDensObj objects. Again, access is through the 
particle ID's. In addition to the individual particles' thermal quantities, the 
total primordial fireball strangeness, baryon, charge, charm, energy, entropy, 
and particle densities, pressure, and Wroblewski factor (see Section 13.7. Ill) are 
included as data members. 

At this level, the constraints on any chemical potentials of constrain-type can 
be imposed, and the correction factors in canonical treatments can be deter- 
mined. Also, as soon as the primordial particle densities are known, the decay 
contributions can be calculated. 



3. 7. 1 Calculating Particle Densities 

Running int GenerateParticleDens clears the current entries in the den- 
sity hash table of the TTMThermalModel object, automatically constrains the 
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chemical potentials (where applicable), calculates the canonical correction fac- 
tors (where applicable), and then populates the density hash table with a 
TTMDensObj object for each particle in the associated set. The decay contribu- 
tions to each stable particle are also calculated, so that the density hash table 
contains both primordial and decay particle density contributions, provided of 
course that the decays have been entered in the associated TTMParticleSet 
object. In addition, the Wroblewski factor and total strangeness, baryon, 
charge, charm and particle densities in the fireball are calculated. 



Note: The summary decay lists of the associated TTMParticleSet object are 
used to calculate the decay contributions. Hence, only stable particles have 
decay contributions reflected in the hash table. Unstable particles that are 
themselves fed by higher-lying resonances, do not receive a decay contribution. 



Each derived class contains the private function int PrimPartDens, which 
calculates only the primordial particle densities and, hence, the canonical 
correction factors, where applicable. In the case of the grand-canonical and 
strangeness-canonical ensembles, this function calculates the densities without 
automatically constraining the chemical potentials of constrain- type first. The 
constraining is handled by int GenerateParticleDens, which calls external 
friend functions, which, in turn, call int PrimPartDens. In the purely canon- 
ical ensemble, int GenerateParticleDens simply calls int PrimPartDens. 
In this way, there is uniformity between the derived classes. Since there is no 
constraining to be done, there is no real need for a separate function in the 
canonical case. 



3.7.2 Calculating Energy and Entropy Densities and Pressure 

GenerateEnergyDens, GenerateEntropyDens and GeneratePressure iterate 
through the existing density hash table and calculate and insert, respectively, 
the primordial energy density, entropy density and pressure of each particle 
in the set. In addition, they calculate the total primordial energy density, 
entropy density and pressure in the fireball, respectively. These functions re- 
quire that the density hash table already be in existence. In other words, int 
GenerateParticleDens must already have been run. If the parameters have 
subsequently changed, then this function must be run yet again to recalculate 
the correction factors or re-constrain the parameters, as required. 
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3.7.3 Bose-Einstein Condensation 



When quantum statistics are taken into account (e.g. in TTMThermalModelBSQ 
or for the non-strange particles in TTMThermalModelBQ), certain choices of pa- 
rameters lead to diverging integrals for the bosons (Bose-Einstein condensa- 
tion). In these classes, a check, based on TTMThermalParticleBSQ : : Parameters- 
Allowed, is included to ensure that the parameters do not lead to problems. 
Including also the possibility of incomplete strangeness and/or charm satu- 
ration (i.e. 75 7^ 1 and/or 7c ^ 1), Bose-Einstein condensation is avoided, 
provided that, 

e (mi-w)/T > 7 ^l 7 ^l ; (38) 



for each boson. If this condition is failed to be met for any of the bosons in 
the set, a warning is issued and the densities are not calculated. 



3.7.4 Accessing the Thermal Densities 

The entries in the density hash table are accessed using the particle ID's. The 
function TTMDensObj* GetDensities(Int_t ID) returns the TTMDensObj ob- 
ject containing the thermal quantities of the particle with the specified ID. The 
primordial particle, energy, and entropy densities, pressure, and decay density 
are extracted from this object using the GetPrimDensity, GetPrimEnergy, 
GetPrimEntropy, GetPrimPressure, and GetDecayDensity functions of the 
TTMDensObj class, respectively. The sum of the primordial and decay particle 
densities is returned by TTMDensObj : : GetFinalDensity. TTMDensObj : :List 
outputs to screen all thermal densities stored in a 
TTMDensObj object. 



ListStableDensities lists the densities (primordial and decay contributions) 
of all those particles considered stable in the particle set associated with the 
model. Access to the total fireball densities is through separate getters defined 
in the TTMThermalModel base class (e.g. GetStrange, GetBaryon etc.). 



3. 7. 5 Further Functions 

GenerateDecayPartDens and GenerateDecayPartDens(Int_t id) (both de- 
fined in the base class) calculate decay contributions to stable particles. The 
former iterates through the density hash table and calculates the decay con- 
tributions to all those particles considered stable in the set. The latter calcu- 
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lates just the contribution to the stable particle with the specified ID. In 
both cases, the primordial densities must be calculated first. In fact, int 
GenerateParticleDens automatically calls GenerateDecayPartDens, so that 
this function does not have to be run separately under ordinary circumstances. 
However, if one is interested in investigating the effect of decays, while keeping 
the parameters (and hence the primordial densities) fixed, then running these 
functions is best (the hash table will not be repeatedly cleared and repopu- 
lated with the same primordial densities). 

ListDecayContributions (Int_t d_id) lists the contributions (in percentage 
and absolute terms) of decays to the daughter with the specified ID. The pri- 
mordial and decay densities must already appear in the density hash table 
(i.e. run int GenerateParticleDens first). ListDecayContribution(Int_t 
p_id, Int_t d_id) lists the contribution of the decay from the specified par- 
ent (with ID p_id) to the specified daughter (with ID d_id) . The percentages 
listed by each of these functions are those of the individual decays to the total 
decay density. 

Next we consider the specific features of the derived TTMThermalModel classes. 
3. 7. 6 TTMTherma IMode IBSQ 

In the grand-canonical ensemble, quantum statistics can be employed and, 
hence, there is a flag specifying whether to use Fermi-Dirac and Bose-Einstein 
statistics or Boltzmann statistics. The constructor, by default, includes both 
the effect of quantum statistics and resonance width. The flags controlling their 
inclusion are set using the SetQStats and SetWidth functions, respectively. 
The functions that calculate the particle, energy, and entropy densities, and 
pressure then use the corresponding functions in the TTMThermalParticleBSQ 
class to calculate these quantities in the required way. The statistics data mem- 
ber (f Stat) of each TTMParticle included in the associated set can be used 
to fine-tune the inclusion of quantum statistics; with the quantum statistics 
flag switched on, Boltzmann statistics are still used for those particles with 
fStat=0. 

In this ensemble, at this stage, both jj, s and /iq can be constrained (either 
separately or simultaneously). In order to accomplish this, the jj, s and/or 
fiQ parameters in the associated TTMParameterSetBSQ object must be set to 
constrain-type. 

It is also possible to constrain ji B by the primordial ratio E/N (the average 
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energy per hadron), n^ + ni (the total primordial baryon plus anti-baryon den- 
sity), or s/T 3 (the primordial, temperature- normalised entropy density). This 
is accomplished by the int ConstrainEoverN, int ConstrainTotalBaryonDens 
and int ConstrainSoverT3 methods, respectively. Running these functions 
will adjust hb such that E/N, n& + ni or s/T 3 , respectively, has the required 
value, regardless of the parameter type of \ib- In addition, the percolation 
model [59] can be imposed to constrain \ib using int ConstrainPercolation. 

This class also accommodates charm, since the associated parameter set in- 
cludes He and ■jc, while the associated particle set may contain charmed 
particles. However, no constraining functions have yet been written for the 
charm content within this ensemble. 

Within the grand-canonical ensemble, it is possible to include excluded volume 
effects. Their inclusion is controlled by the f ExclVolCorrection flag, false 
by default, which is set through the SetExcludedVolume function. When in- 
cluded, these corrections are calculated on calling int GenerateParticleDens, 
based on the hard-sphere radii stored in the TTMParticle objects of the asso- 
ciated particle set. 



3. 7. 7 TTMTherma IMode IBQ 

This class contains the following additional data members: 



f InZtot 


- log of the total partition function, 


flnZO 


- log of the non-strange component of the partition function, 


fExactMuS 


- equivalent strangeness chemical potential, 


fCorrPl 


- canonical correction for S = +1 particles, 


fCorrP2 


- canonical correction for S = +2 particles, 


fCorrP3 


- canonical correction for S = +3 particles, 


fCorrMl 


- canonical correction for S = -1 particles, 


fCorrM2 


- canonical correction for S = -2 particles, 


fCorrM3 


- canonical correction for S = -3 particles. 



Although this ensemble is only applied in the Boltzmann approximation for 
S 7^ hadrons, it is possible to apply quantum statistics to the S = hadrons. 
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This is achieved through the SetNonStrangeQStats function. By default, 
quantum statistics is included for the non-strange hadrons by the construc- 
tors. Resonance width can be included for all hadrons, and is achieved through 
the SetWidth function. The constructors, by default, apply resonance width. 
The functions that calculate the particle, energy, and entropy densities, and 
pressure then use the corresponding functions in the TTMThermalParticle 
classes to calculate these quantities in the required way. 

int GenerateParticleDens populates the density hash table with particle 
densities, including the canonical correction factors, which are also stored in 
the appropriate data members. The equivalent strangeness chemical potential 
is calculated from the canonical correction factor for S — +1 particles. In the 
limit of large VT 3 , this approaches the value of /is in the equivalent grand- 
canonical treatment. 



Running GenerateEntropyDens populates each TTMDensObj object in the 
hash table with only that part of the total entropy that can be unambiguously 
attributed to that particular particle. There is a term in the total entropy that 
cannot be split; this is added to the total entropy at the end, but not included 
in the individual entropies (i.e. summing up the entropy contributions of each 
particle will not give the total entropy). 

At this stage, in this formalism, \xq can be constrained (this is automati- 
cally realised if this parameter is set to constrain- type), while the correlation 
radius (R c ) can be set to the fireball radius (R) by applying the function 
ConserveSGlobally to the associated TTMParameterSetBQ object. 

In exactly the same way as in the grand-canonical ensemble case, /ib can be 
constrained in this ensemble by the primordial ratio E/N (the average energy 
per hadron), n& + ni (the total primordial baryon plus anti-baryon density), 
or s/T 3 (the primordial, temperature- normalised entropy density), as well as 
by the percolation model. 



3. 7. 8 TTMTherrna IModel CanBSQ 

This class contains, amongst others, the following data members: 
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f InZtot - log of the total canonical partition function, 



-PM-.-.D -Pi\yr-.-.o -PM-,-,n 

iMud , iMub , iMuL) 


- equivalent 


chemical potentials, 


f Corrpip 


- correction 


for 


7r^-like Darticles 


f Corrpim 


- correction 


for 


7r~-like particles, 


f Corrkm 


- correction 


for 


i^~-like particles, 


f Corrkp 


- correction 


for 


i^ + -like particles, 


f CorrkO 


- correction 


for 


i^-like particles, 


iUorrakO 


- correction 


for 


i^°-like particles, 


fCorrproton 


- correction 


for 


p-like particles, 


fCorraproton 


- correction 


for 


p-like particles, 


f Corrneutron 


- correction 


for 


n-like particles, 


f Corraneutron 


- correction 


for 


n-like particles, 


f Corrlambda 


- correction 


for 


A-like particles, 


f Corral ambda 


- correction 


for 


A-like particles, 


fCorrsigmap 


- correction 


for 


S + -like particles, 


fCorrasigmap 


- correction 


for 


E _ -like particles, 


iUorrsigmam 


- correction 


for 


X~-like particles, 


f Corrasigmam 


- correction 


for 


S + -like particles, 


fCorrdeltam 


- correction 


for 


A~-like particles, 


fCorradeltam 


- correction 


for 


A~*~-like particles, 


fCorrdeltapp 


- correction 


for 


A ++ -like particles. 


f Corradeltapp 


- correction 


for 


A -like particles. 


f Corrksim 


- correction 


for 


^,~-like particles, 


f Corraksim 


- correction 


for 


H + -like particles, 


fCorrksiO 


- correction 


for 


S°-like particles, 


f CorraksiO 


- correction 


for 


S°-like particles, 


f Corromega 


- correction 


for 


n~-like particles, 


f Corraomega 


- correction 


for 


+ -like particles. 



Since this ensemble is only applied in the Boltzmann approximation, there is 
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no flag for quantum statistics. However, resonance width can be included. This 
is achieved through the SetWidth function. The constructor, by default, ap- 
plies resonance width. The functions that calculate the particle, energy, and 
entropy densities, and pressure then use the corresponding functions in the 
TTMThermalParticle classes to calculate these quantities in the required way. 

int GenerateParticleDens calls int PrimPartDens, which calculates the 
particle densities, including the canonical correction factors, which are then 
also stored in the relevant data members accessible through the double GetCorrFactor 
method. The integrands featuring in the evaluation of the partition function 
and correction factors can be viewed after calling PopulateZHistograms. This 
function populates the array passed as argument with histograms showing 
these integrands as a function of the integration variables (ps and 4>q- Since 
these histograms are created off of the heap, they must be cleaned up after- 
wards. 

GenerateEntropyDens acts in exactly the same way as in the strangeness- 
canonical ensemble case. 



3. 7. 9 Example 

As an example, we consider the strangeness-canonical ensemble, based on the 
particle set and strangeness-canonical parameter set previously defined. After 
instantiating the object, we populate the hash table with primordial and decay 
particle densities: 

root [ ] TTMThermalModelBQ modBQ(&set ,&parBQ) 
root [ ] modBQ . GenerateParticleDens () 
root [ ] parBQ .List () 

***************************** Thermal Parameters **************************** 

Strangeness inside Canonical Volume = 

T = 0.16 (to be FITTED) 

start: 0.16 
range: 0.05 — 0.18 
step: 0.001 

muB = 0.2 (to be FITTED) 

start: 0.2 
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range : — 0.5 
step: 0.001 



muQ 



-0.00636409 



(❖CONSTRAINED*) 



B/2Q: 1.2683 



gammas 



0.8 



(FIXED) 



Can. radius 



6 



(FIXED) 



radius 



6 



(FIXED) 



B/2Q Successfully Constrained 



One notices that the constraint on /iq is now automatically imposed. 

The energy and entropy densities and pressure can be calculated once int 
GenerateParticleDens has been run: 

root [ ] modBQ . GenerateEnergyDens () 
root [ ] modBQ . GenerateEntropyDens () 
root [ ] modBQ . GeneratePressure () 

Now, suppose that we are interested in the thermal densities of the A(1600)° 
and 7T + : 

root [ ] TTMDensObj *delta_dens = modBQ . GetDensities (32114) 
root [ ] delta_dens->List () 

**** Densities for Particle 32114 **** 



n_prim = 0.00138306 
n_decay = 
e_prim = 0.0022912 
s_prim = 0.0139745 
p_prim = 0.000221328 



root [ ] TTMDensObj *piplus_dens = modBQ . GetDensities (211) 
root [ ] piplus_dens->List () 

**** Densities for Particle 211 **** 



n_prim = 0.0488139 
n_decay = 0.119683 
e_prim = 0.0247039 
s_prim = 0.20276 
p_prim = 0.00742708 
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One notices that the 7r + has a decay density contribution, while the A(1600)° 
does not. This is because, unlike the A(1600)°, the n + was considered stable. 



3.7.10 Imposing of Constraints 

The 'Numerical Recipes in C [00] function applying the Broyden globally con- 
vergent secant method of solving nonlinear systems of equations is employed 
by THERMUS to constrain parameters. The input to the Broyden method 
is a vector of functions for which roots are sought. Typically, in the thermal 
model, solutions to the following equations are required (either separately or 
simultaneously): 



(7-1 /t t \ model / 7-> \ collidinq system 

=0 , 
2Q /V J rtr j rncyr( jj a i \2Q J 



primordial 

od 
in 

jji mc "^ e ' / E ^ required 



qmodel qcolliding system r\ 

primordial ' 



^ I ^ J primordial 



0. 



Although, as written, these equations are correct, the quantities B/2Q, S 
and E/N are typically of different orders of magnitude. Since the Broyden 
method in 'Numerical Recipes in C defines just one tolerance level for function 
convergence (TOLF), it is important to 'normalise' each equation: 



(j^ j-y \ model / \ colliding system \ / ^ \ colliding system 

ZQjV) 1 primordial' ' \W J ' W J = 

(omodel _ qcolliding system^ J qcolliding system ri 

\ ^ primordial J ) ' 

' E /V\ mo< ^ e ' / E \ reau ''' re( ^ I / E \ re 1 u ^ rec ^ 



This is the most democratic way of treating the constraints. However, this 
method obviously fails in the event of one of the denominators being zero. For 
the equations considered above, this is only likely in the case of the strangeness 
constraint, where the initial strangeness content is typically zero. In this case, 
where the strangeness carried by the positively strange particles S + is balanced 
by the strangeness carried by the negatively strange particles S_, we write as 
our function to be satisfied, 
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(S/V) 



model 
primordial 



/{\S- 



i I model 
+ \primordial 



/v + \s. 



i I model 
— I primordial 



/v)=o. 



In this way, the constraints can be satisfied to equal relative degrees, and 
equally well fractionally at each point in the parameter space. In addition 
to the constraints listed above, THERMUS also allows for the constraining 
of the total baryon plus anti-baryon density and the temperature-normalised 
entropy density, s/T 3 , as well as the imposing of the percolation model. 



3.7.11 Calculation of the Wroblewski Factor 
The Wroblewski factor [61] is defined as, 



where < uu > + < dd > is the sum of newly-produced uu and dd pairs, while 
all ss pairs are newly-produced if S = in the initial state. 

In THERMUS, A^ is calculated in the following way: 

• Using the primordial particle densities and the strangeness content of each 
particle listed in the particle hash table, the s + s and u + d + u + d densities 
are determined. 

• Assuming S = 0, #s = #s, and so the density of newly-produced ss pairs 
is simply (s + s)/2. 

• From baryon number conservation, the net baryon content in the system, 
tib, originates from the initial state. Thus, 3xn B must correspond to the 
density of u+ d quarks brought in by the colliding nuclei. This is subtracted 
from the total u + d + u + d density to yield the density of newly-produced 
non-strange light quarks. 

• Since #s = #s and, amongst newly-produced non-strange light quarks, 
u + d = u + d, further assuming that fig = implies that u = u = d = d. 
This allows the density of uu and dd pairs to be easily determined. 

3.8 The TTMYield Class 

Often a single experiment releases yields and ratios that contain different 
feed-down corrections. Each yield or ratio then has a different decay chain as- 
sociated with it. Since TTMThermalModel objects allow for just one associated 
particle set, they do not allow sufficient flexibility for performing thermal fits 



As 



2 < ss > 



< uu > + < dd > 



40 



to experimental data. However, TTMThermalFit classes do feature such flex- 
ibility. Before we discuss these classes, let us look at the TTMYield object, 
which forms an essential part of the TTMThermalFit class. 

Information relating to both yields and ratios of yields can be stored in 
TTMYield objects. These objects contain the following data members: 



_C TIT _ 

iName 


- the name of the yield or ratio, 


fIDl 


- the ID of the yield or numerator ID in the case of a ratio, 


fID2 


- denominator ID in the case of a ratio (0 for a yield), 


fFit 


- true if the yield or ratio is to be included in a fit (else predicted), 


fSetl 


- particle set relevant to yield or numerator in case of ratio, 


fSet2 


- particle set relevant to denominator in case of ratio (0 for yield), 


f ExpValue 


- the experimental value, 


f ExpError 


- the experimental error, 


fModelValue 


- the model value, 


f ModelError 


- the model error. 



By default, TTMYield objects are set for inclusion in fits. The functions Fit 
and Predict control the fit-status of a TTMYield object. Particle sets (decay 
chains) are assigned using the SetPartSet method. 

The functions double GetStdDev and double Get QuadDev return the number 
of standard and quadratic deviations between model and experimental values, 
respectively, i.e., 



respectively, while List outputs the contents of a TTMYield object to screen. 
Access to all remaining data members is through the relevant getters and set- 
ters. 



(Model Value — Exp. Value)/Exp. Error, 



(39) 



and, 



(Model Value - Exp. Value)/Model Value, 



(40) 
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3.9 The TTMThermalFit Class 



This is the base class from which the TTMThermalFitBSQ, TTMThermalFitBQ 
and TTMThermalFitCanBSQ classes are derived. Each TTMThermalFit object 
contains: 

• a particle set, the so-called base set, which contains all of the constituents 
of the hadron g well as the default decay chain to be used; 

• a parameter set; 

• a list of TTMYield objects containing yields and/or ratios of interest; 

• data members storing the total x 2 an d quadratic deviation; and 

• a TMinuit fit object. 

A string descriptor is also included in the base class to identify the type of 
model on which the fit is based. This is used, for example, to determine the 
number of parameters in the associated parameter sets. 

Each derived class defines a private function, TTMThermalModel* GenerateThermalModel, 
which creates (off the heap) a thermal model object, based on the base particle 
set and parameter set of the TTMThermalFit object, with the specific quantum 
statistics/resonance width/excluded volume requirements, where applicable. 



3. 9. 1 Populating and Customising the List of Yields of Interest 

The list of yields and/or ratios of interest can be input from file using the 
function InputExpYields, provided that the file has the following format: 



333 


Exp_A 


0.02 




0.01 






-211 


211 


Exp_B 





.990 


0. 


.100 


-211 


211 


Exp_C 





.960 


0. 


.177 


321 


-321 


Exp_C 


1. 


.152 


0. 


.239 



where the first line corresponds to a yield, and has format: 

Yield ID /t Descriptor string /t Exp. Value /t Exp. Error/n 

while the remaining lines correspond to ratios, and have format: 

Numerator ID /t Denominator ID /t Descriptor string /t Exp. Value 
/t Exp. Error/n 
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A TTMYield object is created off the heap for each line in the file, with a 
name derived from the ID's and the descriptor. This name is determined by 
the private function TString GetName, which uses the base particle set to 
convert the particle ID's into particle names and appends the descriptor. In 
addition to all of the PDG ID's in the associated base particle set, the following 
THERMUS-defined identifiers are also allowed: 

• ID = 1: N part , 

• ID = 2: h~, 

• ID = 3: h+. 

A TTMYield object can also be added to the list using AddYield. Such yields 
should, however, have names that are consistent with those added by the 
InputExpYields method; the TString GetName function should be used to 
ensure this consistency. Only yields with unique names can be added to the 
list, since it is this name which allows retrieval of the TTMYield objects from 
the list. If a yield with the same name already exists in the list, a warning is 
issued. The inclusion of descriptors ensures that TTMYield objects can always 
be given unique names. 

RemoveYield(Int_t idl,Int_t id2, TString descr) removes from the list 
and deletes the yield with the name derived from the specified ID's and de- 
scriptor by TString GetName. The TTMYield* GetYield(Int_t idl,Int_t 
id2, TString descr) method returns the required yield. 



3.9.2 Generating Model Values 

Values for each of the yields of interest listed in a TTMThermalFit object arc 
calculated by the function GenerateYields. This method uses the current 
parameter values and assigned particle sets to calculate these model values. 

GenerateYields firstly calculates the primordial particle densities of all con- 
stituents listed in the base particle set. This it does by creating the relevant 
TTMThermalModel object from the base particle set and the parameters, and 
then calling int GenerateParticleDens. In this way, the density hash ta- 
ble of the newly-formed TTMThermalModel object is populated with primor- 
dial densities, as well as decay contributions, according to the base particle 
set (recall that int GenerateParticleDens automatically calculates decay 
contributions in addition to primordial ones). GenerateYields then iterates 
through the list of TTMYield objects, calculating their specific decay contri- 
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butions. New model values are then inserted into these TTMYield objects. In 
addition, the total \ 2 an d quadratic deviation are calculated, based solely on 
the TTMYield objects which are of fit-type. ListYields lists all TTMYield ob- 
jects in the list. 



3.9.3 Performing a Fit 

The FitData(Int_t flag) method initiates a fit to all experimental yields 
or ratios in the TTMYield list which are of fit-type. With f lag=0, a \ 2 fit 
is performed, while f lag=l leads to a quadratic deviation fit. In both cases, 
fit_f unction is called. This function determines which parameters of the 
associated parameter set are to be fit, and performs the required fit using 
the ROOT TMinuit fit class. On completion, the list of TTMYield objects 
contains the model values, while the parameter set reflects the best-fit param- 
eters. Model values are calculated by the GenerateYields method. For each 
TTMYield object in the list, a model value is calculated- even those that have 
been chosen to be excluded from the actual fit. In this way, model predictions 
can be determined at the same time as a fit is performed. ListMinuitlnf o 
lists all information relating to the TMinuit object, following a fit. 

3.9.4 TTMThermalFitBSQ, TTMThermalFitBQ and TTMThermalFitCanBSQ 

The constructor in each of these derived classes instantiates an object with 
the specified base particle set and parameter set and inputs the yields listed 
in the specified file in the TTMYield list. 

The specifics of the fit, i.e. the treatment of quantum statistics (in the grand- 
canonical ensemble and for the non-strange particles in the strangeness-canonical 
ensemble), resonance width (in all three ensembles) and excluded volume cor- 
rections (in the grand-canonical ensemble), are handled through the SetQStats/ 
SetNonStrangeQStats, SetWidth and SetExclVol methods, respectively. By 
default, both resonance width and quantum statistics are included, while ex- 
cluded volume corrections are excluded, where applicable. 



3. 9. 5 Example 

As an example to conclude this section, consider a fit to fictitious particle 
ratios measured in Au+Au collisions at some energy. We will assume a grand- 
canonical ensemble, with the parameters T, /ib and us fitted, and \iq fixed to 
zero. In the grand-canonical ensemble, ratios are independent of the fireball 



44 



radius (this is not true in the canonical ensemble). For this reason, there is no 
need to specify the treatment of the radius. Furthermore, we will ignore the 
effects of resonance width and quantum statistics. 

We begin by instantiating a particle set object, based on the particle list 
distributed with THERMUS. After inputting the particle decays (scaled to 
100%), a parameter set is defined (Note: in the example below $THERMUS must 
be entered in full): 

root [ ] TTMParticleSet set("$THERMUS/particles/PartList_PPB2002.txt") 
root [ ] set. InputDecays("$THERMUS/particles/", true) 
root [] TTMParameterSetBSQ par (0.160, 0.05,0. ,0. ,1.) 

Next, we change the parameters T, \xb and fi$ to fit-type, supplying sensible 
starting values as the arguments to the appropriate functions, as well as the 
range of temperature values for the fit (50 - 180 MeV). For all other properties 
of the fit (step size, fit range etc.), we accept the default values: 

root [ ] par. FitT(0. 160,0.05,0. 180) 
root [ ] par .FitMuB (0.05) 
root [ ] par.FitMuS(0.) 

Next, we prepare a file ('ExpData.txt') containing the experimental data: 



-211 


211 


Exp. 


_A 


0. 


,990 


0. 


.100 


-211 


211 


Exp. 


_B 


0. 


.960 


0. 


,177 


-211 


211 


Exp. 


_D 


1. 


.000 


0, 


,022 


321 


-321 


Exp. 


_B 


1. 


.152 


0. 


.239 


321 


-321 


Exp. 


_D 


1. 


,098 


0, 


,111 


321 


-321 


Exp. 


_C 


1. 


.108 


0, 


.022 


-2212 


2212 


Exp. 


.A 


0. 


,650 


0, 


.092 


-2212 


2212 


Exp. 


_B 


0. 


,679 


0, 


.148 


-2212 


2212 


Exp. 


_D 


0. 


,600 


0, 


.072 


-2212 


2212 


Exp. 


_C 


0. 


,714 


0. 


.050 


-3122 


3122 


Exp. 


_B 


0. 


,734 


0, 


,210 


-3122 


3122 


Exp. 


_C 


0. 


,720 


0, 


.024 


-3312 


3312 


Exp. 


_c 


0. 


.878 


0, 


.054 


-3334 


3334 


Exp. 


_c 


1. 


,062 


0. 


.410 



As one can see, there are multiple occurrences of the same particle-anti- 
particle combination. This is why additional descriptors are required. In this 
case, the descriptors list the particular experiment responsible for the mea- 
surement. In other situations, the descriptors may describe whether feed-down 
corrections have been employed or some other relevant detail that, together 
with the ID's, uniquely identifies each yield or ratio. 
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We are now in a position to create a TTMThermalFitBSQ object based on 
the newly-instantiated parameter and particle sets and the data file. Since 
quantum statistics and resonance width are included by default, we have to 
explicitly turn these settings off: 

root [ ] TTMThermalFitBSQ f it (&set,&par, "ExpData.txt") 
root [ ] fit.SetQStats(kFALSE) 
root [ ] fit.SetWidth(kFALSE) 

Next, let us simply generate the model values corresponding to each of the 
TTMYield objects in the list, based on the current parameters. Part of the 
output of ListYields is shown here: 

root [ ] f it . GenerateYields () 
root [ ] fit. ListYields () 



K+/anti-K+ Exp_B: 



FIT YIELD 

Experiment : 1 . 152 +- 
Model: 0.979833 +- 

Std.Dev.: -0.720365 Quad.Dev 



0.239 


-0.175711 



K+/anti-K+ Exp_D: 



FIT YIELD 
Experiment : 
Model : 
Std.Dev. : - 



1.098 +- 
0.979833 +- 
1.06457 Quad.Dev, 



0.111 


: -0.120599 



K+/anti-K+ Exp_C: 



FIT YIELD 
Experiment 
Model : 
Std.Dev. : - 



1.108 
0.979833 



+- 
+- 



.022 




5.82578 Quad.Dev.: -0.130805 



anti-p/p Exp_A: 



FIT YIELD 
Experiment : 
Model : 
Std.Dev. : - 



0.65 +- 
0.535261 +- 
1.24716 Quad.Dev 



0.092 


: -0.21436 
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Finally, we perform a \ 2 fit: 
root [ ] fit.FitData(O) 



COVARIANCE MATRIX CALCULATED SUCCESSFULLY 

FCN=3. 54326 FROM MIGRAD STATUS=CONVERGED 128 CALLS 129 TOTAL 

EDM=6.43919e-06 STRATEGY= 1 ERROR MATRIX ACCURATE 



EXT PARAMETER 
NO . NAME 

1 T 

2 muB 

3 muS 
EXTERNAL ERROR MATRIX 

4.593e-03 1.289e-03 
1.289e-03 3.689e-04 
5.418e-04 1.548e-04 



VALUE 
1.62878e-01 
3.58908e-02 
1.06828e-02 

NDIM= 
5.418e-04 
1.548e-04 
6.653e-05 



ERROR 
1.10211e-01 
1.91364e-02 
8.13945e-03 
25 NPAR 



STEP 

SIZE 
2.32021e-04 
1.68543e-05 
1.80744e-05 



FIRST 
DERIVATIVE 
-7.98809e-03 
4.05740e-03 
1.25485e-01 



= 3 



ERR DEF=1 



PARAMETER CORRELATION COEFFICIENTS 

NO. GLOBAL 12 3 

1 0.99017 1.000 0.990 0.980 

2 0.99410 0.990 1.000 0.988 

3 0.98814 0.980 0.988 1.000 



FCN=3. 54326 FROM MIGRAD 



STATUS=CONVERGED 



EXT PARAMETER 
NO . NAME 

1 T 

2 muB 

3 muS 



EDM=6.43919e-06 



VALUE 
1.62878e-01 
3.58908e-02 
1.06828e-02 



128 CALLS 129 TOTAL 

STRATEGY= 1 ERROR MATRIX ACCURATE 

PHYSICAL LIMITS 
ERROR NEGATIVE POSITIVE 

1.10211e-01 5.00000e-02 1.80000e-01 
1.91364e-02 . 00000e+00 5.00000e-01 
8.13945e-03 0.00000e+00 5.00000e-01 



Once completed, the associated parameter set contains the best-fit values for 
the fit parameters: 

root [ ] par. List () 

***************************** Thermal Parameters **************************** 



0.162878 



0.110211 



(FITTED ! ) 
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start: 0.16 

range: 0.05 — 0.18 

step: 0.001 



muB 



0.0358908 



0.0191364 



muS 



muQ 
gammas 
radius 

muC 
gammac 



0.0106828 +- 0.00813945 




1 



1 

Parameters unconstrained 



(FITTED ! ) 
start: 0.05 
range : — 0.5 
step: 0.001 

(FITTED ! ) 
start: 
range : — 0.5 
step: 0.001 

(FIXED) 

(FIXED) 

(FIXED) 

(FIXED) 

(FIXED) 



4 Installation of THERMUS 



Having introduced the basic functionality of THERMUS in the previous sec- 
tion, we conclude by outlining the installation procedure. 



Since several functions in THERMUS use 'Numerical Recipes in C code [60] 
(which is under copyright), it is required that THERMUS users have their own 
copies of this software. Then, with ROOT [9] already installed, the following 
steps are to be followed to install THERMUS: 

• Download the THERMUS source; 

• Set an environment variable 'THERMUS' to point at the top-level directory 
containing the THERMUS code; 

• Copy the following 'Numerical Recipes in C [60J functions to $THERMUS/nrc: 
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broydn . 


c 


rsolv . c 


f dj ac . c 


f min . c 


lnsrch. 


c 


nrutil . c 


nrutil . 


h 


qrdcmp . c 


qrupdt . 


c 


rotate . c 


zbrent . 


c; 





• Use the makefiles in $THERMUS/f unctions, $THERMUS/nrc and $THERMUS/ main 
to build the libFunctions . so, libNRCFunctions . so and libTHERMUS . so 

shared object files (run make all in each of these directories); 

• Finally, open a ROOT session, load the libraries and begin: 

root [ ] gSystem->Load(" ./lib/libFunctions.so") ; 
root [] gSystem->Load(" ./lib/libNRCFunctions.so") ; 
root [] gSystem->Load(" ./lib/libTHERMUS.so") ; 
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